﻿using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;
using System.Data;
using System.IO;

namespace PracaInz_v0
{
    class Cluster
    {
        DataTable BasicVect;                     
        DataTable OtherVect;                            // wektory ktore nie weszly do klastra
        public DataRow VectMean;                        // srednia wektorow BasicVect;
        double radius;                                  // promien sredni ze wszystkich odleglosci
        public DataTable SLR;                           // tablica zoptymalizowanych regol logicznych
        int BVectCount;                                 // ilosc wektorow znajdujacych sie w klastrze bazowym
        public int ClusterVectCount;                    // ilosc wektorow w klastrze
        double x = 1.5;                                 // mnoznik "radius" do ustalenia
        public string RegFormula = "Y = ";
        public double NewMean = 0;

        public DataTable AdditionalVect;

        public Cluster(DataTable tab, DataTable rules, string source)
        {
            OtherVect = tab.Copy();
            BasicVect = LoadData.loadFileVect(source);
            AdditionalVect = BasicVect.Clone();

            if (App.ARGS == null)
            {
                SLR = simplyLogicRule(rules);
                BVectCount = BasicVect.Rows.Count;
                mean();
                meanRad();
                distances();
                addVect();
                //scalenie bazowych i dodatkowych wektorow na potrzeby regresji
                BasicVect.Merge(AdditionalVect);
                ClusterVectCount = BasicVect.Rows.Count;
                Regresion();
            }
            else
            {
                BVectCount = BasicVect.Rows.Count;
                mean();
                meanRad();
                distances();
                addVect();
                //scalenie bazowych i dodatkowych wektorow na potrzeby regresji
                BasicVect.Merge(AdditionalVect);
                ClusterVectCount = BasicVect.Rows.Count;
                NewMean = BasicVect.AsEnumerable().Average(n => (double)n.ItemArray[BasicVect.Columns.Count - 1]);
            }

            OtherVect.Dispose();
            BasicVect.Dispose();
        }

        void mean()                 // liczy srednia
        {
            VectMean = BasicVect.NewRow();

            for (int i = 0; i < BasicVect.Columns.Count; i++)           // ostatnia to "y"
            {
                double tmp = 0;
                for (int j = 0; j < BasicVect.Rows.Count; j++)
                {
                    tmp += Convert.ToDouble(BasicVect.Rows[j][i]);
                }
                VectMean[i] = (tmp / BasicVect.Rows.Count);
            }
        }

        void meanRad()                  // liczy srednia z odleglosc od srodka do wektorow bazowych w klastrze
        {
            double tmp = 0;

            double dist = 0;

            for (int i = 0; i < BasicVect.Rows.Count; i++)                      // liczy odleglosci od srodka klastra
            {
                for (int j = 0; j < BasicVect.Columns.Count - 1; j++)           // ostatnia to "y"
                {
                    tmp += Math.Pow((Convert.ToDouble(BasicVect.Rows[i][j]) - Convert.ToDouble(VectMean[j])), 2);
                }
                tmp = Math.Sqrt(tmp);

                dist += tmp;

                tmp = 0;
            }

            radius = (dist / BasicVect.Rows.Count) * x;
        }

        void distances()                                                        //oblicza odleglosci
        {
            IEqualityComparer<DataRow> comparer = DataRowComparer.Default;

            for (int i = 0; i < BasicVect.Rows.Count; i++)                      // porownuje i usuwa niepotrzebne wektory
            {
                for (int j = 0; j < OtherVect.Rows.Count; j++)
                {
                    bool bEqual = comparer.Equals(BasicVect.Rows[i], OtherVect.Rows[j]);

                    if (bEqual)
                    {
                        OtherVect.Rows[j].Delete();
                    }
                }
            }

            OtherVect.Columns.Add("dist", typeof(double));                      // dodaje kolumne z odleglosciami od srodka

            double tmp = 0;

            for (int i = 0; i < OtherVect.Rows.Count; i++)                      // liczy odleglosci od srodka klastra
            {
                for (int j = 0; j < OtherVect.Columns.Count - 2; j++)           // -2 bo ostatnia to "dist" a przedostatnia "y" 
                {
                    tmp += Math.Pow((Convert.ToDouble(OtherVect.Rows[i][j]) - Convert.ToDouble(VectMean[j])), 2);
                }
                tmp = Math.Sqrt(tmp);

                OtherVect.Rows[i]["dist"] = tmp;

                tmp = 0;
            }

            DataView dv = OtherVect.DefaultView;
            dv.Sort = "dist ASC";                                               // sortuje tabele od najblizszych do najdalszych
            OtherVect = dv.ToTable();
        }

        public DataTable simplyLogicRule(DataTable rules)                       // upraszcza regoly logiczne
        {
            DataTable SLR = rules.Clone();                                      // tablica uproszczonych regol
            DataTable copy = rules.Copy();                                      // kopia 

            for (int i = 0; i < copy.Rows.Count; i++)
            {
                double minVal = Convert.ToDouble(copy.Rows[i][2]);              // jeśli znak '<' szukaj najmniejszej wartosci
                double maxVal = Convert.ToDouble(copy.Rows[i][2]);              // jeśli znak '>' szukaj najwiekszej wartosci

                for (int j = i + 1; j < copy.Rows.Count; j++)
                {
                    if (copy.Rows[i][0].Equals(copy.Rows[j][0]) && copy.Rows[i][1].Equals(copy.Rows[j][1]))
                    {
                        if (Convert.ToChar(copy.Rows[i][1]) == '<')
                        {
                            if (minVal > Convert.ToDouble(copy.Rows[j][2]))
                                minVal = Convert.ToDouble(copy.Rows[j][2]);
                        }
                        else        // znak == '>'                        
                        {
                            if (maxVal < Convert.ToDouble(copy.Rows[j][2]))
                                maxVal = Convert.ToDouble(copy.Rows[j][2]);
                        }
                        copy.Rows[j].Delete();
                        j--;                                                    // gdy usuwamy musimy sie cofnac
                    }
                }

                if (Convert.ToChar(copy.Rows[i][1]) == '<')
                {
                    copy.Rows[i][2] = minVal;
                }
                else        // znak == '>'
                {
                    copy.Rows[i][2] = maxVal;
                }

                SLR.Rows.Add(copy.Rows[i].ItemArray);
            }

            return SLR;
        }

        void addVect()                                                  // dopisuje najblizsze wektory z othervect do basicvect 
        {
            double r = 0;                                               // odleglosc od srodka klastra dopisywanego wektora
            DataRow dr = BasicVect.NewRow();

            while (r < radius && AdditionalVect.Rows.Count < BVectCount )
            {
                for (int i = 0; i < BasicVect.Columns.Count; i++)
                {
                    dr[i] = OtherVect.Rows[0][i];
                }
                r = Convert.ToDouble(OtherVect.Rows[0][OtherVect.Columns.Count - 1]);
                OtherVect.Rows[0].Delete();
                AdditionalVect.Rows.Add(dr.ItemArray);
            }
        }

        void Regresion()
        {
            if (BasicVect.Rows.Count > 2)
            {
                double[,] DX = new double[BasicVect.Rows.Count, BasicVect.Columns.Count - 1];

                for (int i = 0; i < BasicVect.Rows.Count; i++)
                {
                    for (int j = 0; j < BasicVect.Columns.Count - 1; j++)
                    {
                        DX[i, j] = Convert.ToDouble(BasicVect.Rows[i][j]);
                    }
                }

                double[,] DY = new double[BasicVect.Rows.Count, 1];

                for (int i = 0; i < BasicVect.Rows.Count; i++)
                {
                    DY[i, 0] = Convert.ToDouble(BasicVect.Rows[i]["Y"]);
                }

                //// B = (X'X)^-1X'y
                double[,] XT = Transpose(DX);
                double[,] XTX = MultiplyMatrix(XT, DX);
                double[,] XI = Invert(XTX);
                double[,] XY = MultiplyMatrix(XT, DY);
                double[,] B = MultiplyMatrix(XI, XY);                   // tablica wspolczynnikow rownania

                for (int i = 0; i < B.GetLength(0); i++)                // string z rownaniem
                {
                    RegFormula += Math.Round(B[i, 0], 5).ToString() + " * x" + (i + 1).ToString();
                    if (i < B.GetLength(0) - 1)
                    {
                        RegFormula += " + ";
                    }
                }                

                NewMean = BasicVect.AsEnumerable().Average(n => (double)n.ItemArray[BasicVect.Columns.Count - 1]);
            }
            else
            {
                NewMean = BasicVect.AsEnumerable().Average(n => (double)n.ItemArray[BasicVect.Columns.Count - 1]);
            }

        }

        private double[,] Transpose(double[,] a)                        // transponowanie macierzy
        {
            double[,] t = new double[a.GetLength(1), a.GetLength(0)];

            for (int i = 0; i < a.GetLength(0); i++)
            {
                for (int j = 0; j < a.GetLength(1); j++)
                {
                    t[j, i] = a[i, j];
                }
            }

            return t;
        }

        private double[,] MultiplyMatrix(double[,] a, double[,] b)        // mnozenie macierzy
        {
            double[,] c = new double[a.GetLength(0), b.GetLength(1)];
            if (a.GetLength(1) == b.GetLength(0))
            {
                for (int i = 0; i < c.GetLength(0); i++)
                {
                    for (int j = 0; j < c.GetLength(1); j++)
                    {
                        c[i, j] = 0;
                        for (int k = 0; k < a.GetLength(1); k++) // OR k<b.GetLength(0)
                            c[i, j] = c[i, j] + a[i, k] * b[k, j];
                    }
                }
            }
            return c;
        }


        /// <summary>
        /// Calculate the inverse of a matrix using Gauss-Jordan elimination.
        /// </summary>
        /// <param name="data">The matrix to invert.</param>
        /// <param name="inverse">The inverse of the matrix.</param>
        /// <returns>false of the matrix is singular, true otherwise.</returns>
        private double[,] Invert(double[,] data)
        {
            double[,] inverse = new double[data.GetLength(0), data.GetLength(1)];

            if (data == null)
                throw new ArgumentNullException("data");
            if (inverse == null)
                throw new ArgumentNullException("null");

            int drows = data.GetLength(0),
                dcols = data.GetLength(1),
                irows = inverse.GetLength(0),
                icols = inverse.GetLength(1),
                n, r;
            double scale;

            //Validate the matrix sizes
            if (drows != dcols)
                throw new ArgumentException(
                "data is not a square matrix", "data");
            if (irows != icols)
                throw new ArgumentException(
                "inverse is not a square matrix", "inverse");
            if (drows != irows)
                throw new ArgumentException(
                "data and inverse must have identical sizes.");

            n = drows;

            //Initialize the inverse to the identity
            for (r = 0; r < n; ++r)
                for (int c = 0; c < n; ++c)
                    inverse[r, c] = (r == c) ? 1.0 : 0.0;

            //Process the matrix one column at a time
            for (int c = 0; c < n; ++c)
            {
                //Scale the current row to start with 1

                //Swap rows if the current value is too close to 0.0

                if (Math.Abs(data[c, c]) <= 2.0 * double.Epsilon)
                {
                    for (r = c + 1; r < n; ++r)
                        if (Math.Abs(data[r, c]) <= 2.0 * double.Epsilon)
                        {
                            RowSwap(data, n, c, r);
                            RowSwap(inverse, n, c, r);
                            break;
                        }
                    //if (r >= n)
                    //    return false;
                }
                scale = 1.0 / data[c, c];
                RowScale(data, n, scale, c);
                RowScale(inverse, n, scale, c);

                //Zero out the rest of the column

                for (r = 0; r < n; ++r)
                {
                    if (r != c)
                    {
                        scale = -data[r, c];
                        RowScaleAdd(data, n, scale, c, r);
                        RowScaleAdd(inverse, n, scale, c, r);
                    }
                }
            }

            return inverse;
        }

        /// <summary>
        /// Swap 2 rows in a matrix.
        /// </summary>
        /// <param name="data">The matrix to operate on.</param>
        /// <param name="cols">
        /// The number of columns in <paramref name="data"/>.
        /// </param>
        /// <param name="r0">One of the 2 rows to swap.</param>
        /// <param name="r1">One of the 2 rows to swap.</param>
        private void RowSwap(double[,] data, int cols, int r0, int r1)
        {
            double tmp;

            for (int i = 0; i < cols; ++i)
            {
                tmp = data[r0, i];
                data[r0, i] = data[r1, i];
                data[r1, i] = tmp;
            }
        }

        /// <summary>
        /// Perform scale and add a row in a matrix to another
        /// row:  data[r1,] = a*data[r0,] + data[r1,].
        /// </summary>
        /// <param name="data">The matrix to operate on.</param>
        /// <param name="cols">
        /// The number of columns in <paramref name="data"/>.
        /// </param>
        /// <param name="a">
        /// The scale factor to apply to row <paramref name="r0"/>.
        /// </param>
        /// <param name="r0">The row to scale.</param>
        /// <param name="r1">The row to add and store to.</param>
        private void RowScaleAdd(double[,] data, int cols, double a, int r0, int r1)
        {
            for (int i = 0; i < cols; ++i)
                data[r1, i] += a * data[r0, i];
        }

        /// <summary>
        /// Scale a row in a matrix by a constant factor.
        /// </summary>
        /// <param name="data">The matrix to operate on.</param>
        /// <param name="cols">The number of columns in the matrix.</param>
        /// <param name="a">
        /// The factor to scale row <paramref name="r"/> by.
        /// </param>
        /// <param name="r">The row to scale.</param>
        private void RowScale(double[,] data, int cols, double a, int r)
        {
            for (int i = 0; i < cols; ++i)
                data[r, i] *= a;
        }
    }
}
